#install.packages("maptools")
#install.packages("maps")
#install.packages("rgdal")
#install.packages("rgeos")
#install.packages("ggmap")
#install.packages("ggplot2")
#install.packages("plyr")
#install.packages("grid")
#install.packages("gridExtra")
library(plyr)
library(ggplot2)
library(ggmap)
library(maps)
library(maptools)
if (!require(gpclib)) install.packages("gpclib", type="source")
gpclibPermit()
library(rgdal)
library(rgeos)
library(grid)
library(gridExtra)
block_url<-"http://www.mdp.state.md.us/msdc/census/cen2010/maps/tiger10/blk2010.zip"
download.file(block_url,destfile="./data/blk2010.zip")
dir.create("./data/blk2010")
unzip("./data/blk2010.zip",exdir="./data/blk2010")
map <- readOGR("./data/blk2010", "blk2010")
map_wgs84 <- spTransform(map, CRS("+proj=longlat +datum=WGS84"))
bal.blk<-map_wgs84[map_wgs84@data$CNTY2010==24510,]
save(bal.blk,file="./data/bal_block.rda")
neighbor_url<-"http://bniajfi.org/wp-content/uploads/2016/04/VS-14-Health.zip"
download.file(neighbor_url,destfile="./data/neighbor.zip")
unzip("./data/neighbor.zip",exdir="./data/neighbor")
map.neighbor<-readOGR("./data/neighbor/VS 14 Children & Family Health","VS14_Health")
bal.nei <- spTransform(map.neighbor, CRS("+proj=longlat +datum=WGS84"))
save(bal.nei,file="./data/bal_neighbor.rda")
load("./data/bal_block.rda")
load("./data/bal_neighbor.rda")
CenterOfMap <- geocode(" 39.299768,-76.614929")
Baltimore <- get_map(c(lon=CenterOfMap$lon, lat=CenterOfMap$lat),zoom = 12, maptype = "road", source = "google")
BaltimoreMap <- ggmap(Baltimore)
BalHeatmap<-function(shpdata){
shpdata@data$id = rownames(shpdata@data)
shpdata.points = fortify(shpdata, region="id")
shpdata.df = join(shpdata.points, shpdata@data, by="id")
bal_heatmap<- BaltimoreMap+
geom_polygon(aes(long,lat,group=group,fill=LifeExp14),data=shpdata.df,alpha=0.8,color="#999999",size=0.2)+
scale_fill_gradient(low = "#ffffcc", high = "#ff4444", space = "Lab", na.value = "grey50",guide = "colourbar")
return(bal_heatmap)
}
# imput neihgborhood LE14 to block level
bal.blk@data$LifeExp14<-over(bal.blk,bal.nei)$LifeExp14
#BalHeatmap(bal.blk) # plot LE4 on block level
# creat data frame (tidy data)
smooth_bal<-function(dataline,dataframe,lim){
subdata<-dataframe[ dataframe[,2]>=(as.numeric(dataline[2])-lim) &
dataframe[,2]<=(as.numeric(dataline[2])+lim) &
dataframe[,3]>=(as.numeric(dataline[3])-lim)
& dataframe[,3]<=(as.numeric(dataline[3])+lim) , ]
new_var<-mean(subdata[,4],na.rm=TRUE)
return(new_var)
}
bal.blk@data$LifeExp14<-over(bal.blk,bal.nei)$LifeExp14
data.blk<-data.frame(blk=bal.blk@data$BLK2010,long=coordinates(bal.blk)[,1],lat=coordinates(bal.blk)[,2],LifeExp14=bal.blk@data$LifeExp14)
bal.blk@data$LifeExp14<-apply(data.blk,1,FUN=smooth_bal,dataframe=data.blk,lim=0.015)
Map.nei<-BalHeatmap(bal.nei)+
ggtitle("Fig2: Baltimore City Life Expectancy on Neighborhood Level")+
theme(text = element_text(size=12,face="bold"))
# LE4 on neighborhood level
Map.blk<-BalHeatmap(bal.blk)+
ggtitle("Fig3: Baltimore City Life Expectancy on Block Level")+
theme(text = element_text(size=12,face="bold"))
# plot LE4 after smoothing on block level
The inequality in Baltimore has been a big issue: “Roland Park would be the 4th longest-living country in the world, while Seton Hill would be the 230th. Fourteen Baltimore neighborhoods have lower life expectancies than North Korea. Eight are doing worse than Syria.” (Christopher Ingraham, The Washington Post, 2015). In this study, we focus on life expectancy on block level in Baltimore City. Our goal is to develop a model for predicting life expectancy in Baltimore down to single block resolution with estimates of uncertainty.
For geographic information about blocks in Baltimore, we use the Tiger shape files from Census 2010, from which we obtained block ID, coordinates and footprint for each block. There are 13488 blocks defined by Census 2010, and we use the polygon center of each block to calculate the map distance between two locations (use Google Maps API through R package “ggmap”).
#fellspoint<-bal.nei[bal.nei@data$CSA2010=="Fells Point",]
fpblk<-bal.blk[bal.blk@data$BG2010%in%c(245100203001,245100203002,245100203003,245100202001,245100202002,245100201001,245100201002,245100105001,245100105002),]
plot(fpblk,axes=TRUE,main="Fig1: Block Boundaries and centers of Fells Point Neighborhood")
points(coordinates(fpblk),pch=1,cex=0.7)
legend("bottomright",legend=c("Block center","block boundary"),cex=0.7,lwd=c(NA,1),pch=c(1,NA),text.font = 10)
We could only obtain life expectancy data for Baltimore City on neighborhood level (http://bniajfi.org/vital_signs/data_downloads/) through available resources. Therefore, we need to generate reliable outcome variable for 13488 blocks according to the life expectancy data for 55 neighborhoods in Baltimore City. We applied a 2 dimensional smoothing method to generate the outcome variable for Baltimore City’s blocks.
grid.arrange(Map.nei,Map.blk,nrow=1)